The fractional Brownian motion property of the 
turbulent refractive index and the Fermat's 
Extremal Principle 

Dario G Perez|J 

Centre de Investigaciones Opticas (CIOp), CC. 124 Correo Central, 1900 La Plata, 
Argentina. 

E-mail: dariop@minet.uni-jena.de 

Abstract. Wc introduce fractional Brownian motion processes (fBm) as an 
alternative model for the turbulent index of refraction. These processes allow to 
reconstruct most of the index properties, but they are not differentiable. We overcome 
the apparent impossibility of their use in the variational equation coming from 
the Fermat's Principle with the introduction of a Stochastic Calculus. Afterwards, 
we successfully provide a solution for the stochastic ray equation; moreover, its 
implications in the statistical analysis of experimental data is discussed. 

PACS numbers: 02.50.Fz, 42.15.-i, 42.25.Dd, 47.27.Eq. 
1. Introduction 

Diverse experimental techniques have been devoted to the study of the optical properties 
of the turbulent atmosphere. Most of these techniques are based on the analysis of the 
ouput of laser beams making their way through it. But also, controlled experiences had 
been developed for the laboratory, such as the experiments performed by Consortini et 
al 0, D, These experiences apply Geometric Optics to interpret the data acquired. 
This analysis has its theoretical grounds on the precursor paper by P. Beckman who 
was able to find a nice relationship between the variance of the turbulent refractive index 
/.i(r)-being homogeneous and isotropic-and the variance of the laser beam wandering 
over a screen. This derivation, however, is based upon the strong assumption that the 
stochastic process is smooth enough to allow continuity in its derivatives. 

This assumption has serious problems when checked against Tatarski's foundational 
work about lightwave propagation through turbulent atmosphere [|5|. A short revision 
of the covariance of the turbulent index of refraction found by him shows that the 
stochastic process associated is nowhere differentiable p. 

However, the model introduced by Beckmann is just one of the many used during the 
last 40 years to solve imaging problems through the atmosphere, with relative success. 

f This work was carried out while visiting the Mathematisches Institut at the Friedrich-Schiller- 
Universitat in Jena, Germany. 
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In fact, it is clear that the stochastic properties of the refractive index in Atmospheric 
Optics have never been fully understood nor explained. For instance, most works treated 
it like a Gaussian process ^ , some others suggested stationary increments, while other 
works have proposed the use of non-Gaussian statistics. Ishimaru's book p presents an 
extensive description of these works. 

Simultaneously, over the last decade an intense debate in Fluid Dynamics has been 
carried out about whether or not passive scalar fields, among which is the turbulent 
refractive index, behave like the velocity field. That is, under which circumstances they 
inherit the stochastic properties of the turbulent velocity u(r), which is a Gaussian 
process that follows the Kolmogorov refined similarity hypotheses inside the inertial 
range, 

E[||Au(r)r] = Ar^E[e^/^] Hrf/^, k « ||r|| « Lq, (1) 

where Er is the dissipated energy. An is a constant depending on n. Iq {inner length) and 
Lq [outer length) are constants that determine the inertial range, and can be estimated 
theoretically. When intermittence effects are noticeable the dissipated energy modifies 
slightly the right-hand side of this equation, changing the power over the correlation 
distance 

(2) 

here Cn is called multi-fractal exponent. 

It was shown |T^, |TT| that passive scalars fields are nearly gaussian as far as Iq <^ Lq, 
and are unsatisfactory modeled either by log-normal distributions, or the so called 
Frisch's /9-model. Effectively, in the inertial range, these fields obey a law resembling 
the Kolmogorov's law for the velocity fields. Moreover, they do not present the multi- 



fractal property due to intermittence whenever an isotropic velocity field is present [12 



i.e. Cn = 0. Hence, the behavior of the turbulent refractive index predicted by Tatarski 
has been confirmed and extended. 



In the meantime, Stolovitzky and Sreenivasan |jT3[ successfully obtained 
equation (|l|) modeling the turbulent velocity field as a fractional Brownian motion 
(fBm). But this model failed to replicate the intermittent property. It must be stressed 
that Kolmogorov refined similarity hypotheses also implies that the velocity field is 
independent of the dissipated energy probability distribution. 

According to what we pointed out here, the fBm processes seem to be a good 
alternative model for the turbulent refractive index. First, they are gaussian, and second, 
they let us test the Structure Function's power factor. Also, they are continuos but 
nowhere differentiable, as it results from the application of the Kolmogorov hypotheses 
to the refractive index. In few words, the fractional Brownian motion model for 
the turbulent index of refraction describes closer the turbulent properties of passive 
scalar fields than Beckmann's proposed model. However, a mathematical difficulty is 
introduced since we have lost differentiability in the usual sense, and as a consequence 
we cannot appeal to the variational methods used in Geometric Optics. Our task here 
will be to show how using Stochastic Calculus techniques (for a first introduction we 
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refer the reader to 0ksendal |T^, 0) this situation can be overridden providing an 
exphcit solution to the ray propagation through air in turbulent motion. Besides the 
intrinsic complexity of these tools, our model is meant to provide a bridge between the 
stochastic processes and the experimental data. Also, we will gain knowledge about 
problems handling (geodesic) stochastic equations in Optics. 

2. Stochastic Differential Equations in Geometric Optics 

2.1. The Ray- Path Equations 

Fermat's Extremal Principle in Geometric Optics, means to find the variational solution 
to 

6(^j nds^ = 0. (3) 

The solution to this equation is interpreted as the ray trajectory, which we shall denote 
here as q{T). The parameter r has, in principle, no physical meaning. In Optics Treatises 
this parameter is usually replaced with one of the trajectory coordinates, which fulfills 
dqi/dr > 0, and is thus called the propagation direction. But the election of this 



parameter can not be done at will [|T^], since, for any parameterization chosen, the 
Optical Lagrangian 

L{q,q) =n{q)\\q\\, (4) 

(g, g G are the position and velocity respectively) is degenerated. Its solution is not 
univocally determinated because 



for any pair {q,q). That is, calculating the momentum, 
we can write the Lagragian as follows, 

^(^'^) = E§^^ = E^^^- (7) 

i i 

This is homogeneous in the velocities which allows us to recalculate the momentum and 
find, 

Qq] z—j y Qqi dq'dq^ J 

"-E^^'. Vf. (8) 

SO this matrix is singular as we stated above. Another consecuence can be derived from 
equation (j^); also, it induces the following relation 

\\pf = n\q), (9) 
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which indicates that the choice of coordinates and momenta is not free. 

Now, the degeneracy of the Lagragian can be worked out in the Hamiltonian 
framework using the constraint we have just found. This problem of constrained 
Hamiltonians is known as Dirac's problem in the literature The procedure is 

as follows: first, define the constraint function 

nP,Q) = l[\\pf-n\q)]; (10) 

then, calculate the Hamiltonian H from the original Lagragian, from equations (|]) and 
(0) it is exactly zero 

H = -L = J2p^<i' - T.P^'i' = 05 (11) 

i i i 

finally, build a new Hamiltonian 

H{p, q) := H{p, q) + A^(p, q) = A^(p, q) (12) 

and apply the variational procedure to all the coordinates included A. 
By doing so, we obtain the following dynamic equations 

dH , dm A. 



(13) 



^ dq ^ dq 2^'^ 

dH am ^ 

op op 
and the constraint, over the phase space, 

= vl/(p,g) = i[||pf -n2(g)]. (14) 

To ensure A is well defined we need to introduce compatibility conditions. Which arise 
from establishing the constraints as motion constants, that is, / = = {/, H} with {■, ■} 
the Poisson bracket. In our problem these conditions reduce just to one: {H, = 
which is automatically accomplished. There are no secondary constraints derivated from 
the compatibility conditions so ( pTsD and ( p^lf ) completely define our problem . Notice 
that A is actually a smooth function on the constrained space that can be freely chosen. 
Finally, this pair of equations can be combined to yield 

dr \\dTj 2 ^ (15) 
||gf = AV(g). 

We now realize that with each selection we make for A the parameter r is also set, i.e. 
if we choose A = then 

= 1 and ds = \\q\\dT = dr (16) 

r is then the arc-length. But selecting A = 1 gives us ds = ndr and now the parameter 
is r = J ds/n. 
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2.2. Linearizing the Trajectory Equations 



The ray equations we have just found are nonhnear, so in this section we are going to 
hnearize them and at the same time define the parameter r. Therefore, let n be the 
refractive index of the medium and uq its average, we write 



■Tin 



«e(g). 



(17) 



e(g) represents a perturbation field with its intensity measured by a. We also assume 
that it contains all the inhomogeneities of the media, so when a = the index is 
constant. Now we express the solution to (|l^) in a power series on a: 

qir) = qo + aqi + a'^q2^ . (18) 

Although we can also develop a series for the constraint function A, it is far more 
convenient to set its value beforehand. Let us rewrite the first equation in ([TSl) as follows 

From all the possible parametrizations we choose (|16|), so 

1 _ 1 ^ (-l)"6"(g) . 



n'^{q) 



n=l 



n. 



2n+2 




a 



(19) 



(20) 



in short we will write := (— l)"e"(g)/nQ'^'''^. 

Inserting the former series for A^ in (|l^), after some algebraic manipulation we 
obtain the following family of differential equations 
d'^qo 



dr^ 
(fqi 

dr^ 
dr^ 



0, 




' dr J dr 



and when n > 1 : 



(21) 



E 

m=l 



m / 



dq. 



m~k 



dr 



dq„ 



dr 



the constraint condition also gives us a constraint equation for each differential equation 
in(0); 

2 

= 1, 



dqo 
dr 
dqi dqo 
dr dr 

E dqk dqn-k 
dr dr 

k=0 



0, 



(22) 



0, 



for all n > 1. 



We can readily find the zero-order solution of the first equation in (^). The result 
is the linear relationship: qo{r) = ar + b. Given that the boundary condition to this 
problem is 

g(0) = 0, (23) 
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it implies that b = 0. Now we use the constraint condition to obtain ||ap = 1, so we 
are free to choose the coordinate frame best suited to our purposes. Let this coordinate 
frame be: 



qo = re„ (24) 

this will be our forward direction of propagation. To solve the remaining equations we 
need a bit more than algebra. 

The turbulent refractive index measures the separation between the index of 
refraction and its average; /i(r) := n(r) — no. It is a small quantity, that is how its 
increments are often replaced in the literature by those of the permitivity. This passive 
scalar field follows the Kolmovorv refined hypotheses in the inertial range, and so its 
Sturcture Function is 



Deir) 



E 



(n^(r + 



72^1 r 



'))' 



4:Ctl. 



2,2/3-2//,, II 2// 



(25) 



is the permitivity structure constant and H is some positive constant less than one. 
If the turbulence is isotropic and homogeneous then the Kolmogorov hypotheses sets 
H = 1/3; so, we have introduced the inner scale to correct the departure from this 
ideal situation. Then, according to what we have just said at the introduction and the 
definition (0), we propose the following: 

e(r) ■.= B^l,^r\\) 

is a fractional Brownian motion, and H G (0, 1) is the Hurst parameter. 
Gaussian process with the following properties [pO|] : 



(26) 
It is a 



E[B 



Hi 



1 



E[B^{s)B^{t)\ = - 
and the self- similarity property 



\2H 



0, (27) 
+ \tf" - \s - tf"] , (28) 

(29) 

this last equation implies that both variables have the same probability distribution. 
Thus from equations (|25| ) and (^6]) we have. 



B 



as 



HuHl 



a"B 



2,2/3-2// |, 1 1 2// 



£'■0 



a^l. 



{B^il 



H li-l I 
I 



\2H 



whenever ||r|| ^ 11 r' 



■'0 

so it is, 



a 



21 



1/3 



'CI 



(30) 



(31) 



Estimates for the Structure Constant and the inner length tell us that a ~ 10^^. 
Therefore, in order to examine the stochastic behavior of a wandering beam will be 
enough to consider the first order solution. 

The first order constraint condition reads then 
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and together with boundary condition (^) makes the component along the 2;-axis null 
all over the ray trajectory. Of course, this condition agrees with the corresponding 
dynamical equation ( pT)) . So we are left with a differential equation for the perpendicular 
displacements to the direction of propagation. Finally, we multiply these displacements 
by a, and write the first-order equation as 



^Q = ^V,e(ze. + Q). 



(33) 



We must provide a context to understand the previous equation. That is, a 
stochastic equation is not only determined by the type of process (the fractional 
Brownian motion in our case) attached to it, but also by the integro- differential theory 
employed to define its derivatives. Moreover, there are distinctive stochastic integration 
methods whether H > 1/2 oi H < 1/2 ||21|. Here we are going to make use of 
the Stochastic Calculus exposed in the appendix, so only the H > 1/2 case will be 
considered. By doing so, either we are considering the inertial-diffusive range, in the 
following sense 



1 

U> 3, 



(34) 



or the anisotropic scalar situation 1 |22|. This situation could be observed in a 



laboratory if an isotropic velocity field can not achieved by the experimental setup. 

Because the turbulent refractive index oscillates around its mean value, it is 
expected that the light wanders around the z-axis over the screen. So the solution we 
are looking for must have expectation zero. This can easily achieved by the formalism 
we are employing: the stochastic integrals (formally known as fractional ltd integrals) 
defined by the fractional white noise and Wick product on fractional Hida spaces have 
expectation zero. Henceforth, from definition (^) we can calculate the gradient of the 
index of refraction; we have, using the continuity and differential properties for the 
fractional Brownian motion-equations ( |A.27|) and ( |A.3(]| )-in and applying the chain 
rule, the following 



d 



dB 



H 



ds 



X 



_ W^^(/o^||r||) 



la 



X 



1,2,(35) 



is the fractional white noise. Remember, once more, this last identity must be 
understood in terms of the formal definition of white noise inside the fractional Hida 
spaces, it has nothing to do with the usual concept of derivative. 

Next, the procedure to interpret equation (|33|) requieres to replace all the ordinary 
products containing stochastic variables by Wick products, and we finally write 

'W{lo'\\ze, + Q\\)' 



rf2 



a 



2/0^0 



ze, 



+ Q\\ 



oQ. 



(36) 



The fractional white noise is a functional on the real line and its composition with 
another stochastic process has to be defined. Because any analytic function is expressed 
by a power series, as 0ksendal et al [O suggest, we follow our sustitution rule for 
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products and relace the powers in the series by Wick's powers-just as we did with the 
Wick exponential in the appendix ( |A. 33| )-whenever a stochastic process is an argument 
for the given function. The representation for the noise in S"^ is a series with analytic 
functions as components ( [A.29| ); thus, we change these components 



z) Us) ds ^ / <P%s, Z) Us) ds, (37) 
Jr 

Z is some continuous stochastic process with K[Z] := zq ^ 0, and 0*(s, ■) is the Wick 
representation of 0(s, ■) = H{2H — 1) |s — - ^^'"^ . We can write H^re^ + <5|| — z + Q'^ /2z 
because Q ~ 0{a), and then we can evaluate the fractional white noise aX z + a^Ziuj): 



[s, z + a^Z) = H{2H - 1) U + a^Z 



]o{2H-2) 
o{2H-2) 



= H{2H - 1) [{z - s) + a'zX''''-'' , (38) 

we have just took the positive part of the absolute value; it is enough for us examine 
this situation. If E[a;^Z'] = a'^ZQ then 

[{z-s) + a'zY^'^''^=[{z-s) + a'zor-' + 

, ^ a^-i2H-2)...{2H-3-n) (39) 
n! [(.-.) + «%]— ' 

and all the terms in the series are of order higher or equal to 2 in a. We just need to 
compare the first term against the deterministic coeficient in the white noise series, 

(j){s, t + o?zq) - (f){s, t) ~ (t - s)^^^-^^a^ (40) 

and because this happens coordinate to coordinate in the fractional white noise 
decomposition we find 

^-•"('0'^) _ ^-^fe'll'^'+QII) . O(a') (41) 
z \\~e, + Q\\ ^ ^ ' 

The first-order equation (^) is unaffected by this replacement since they differ in a^. 
We have arrived to the linear equation 

^ r^f.^ _ . W''{l^H) o Q{z) 
dz 

we have set g = q;/2/o^o- 

3. The Stochastic Volterra Equation 

3.1. The Stochastic Volterra Equation and Its Solution 
The integral form of equation ( ^2]) is, 

(43) 



'■5 

^Q{z)=g ■ -/O-vv-. ^ ^42) 



^0 



Qiz)=Qoz + g I I '-^ '-oQ{s')ds'ds 



Qoz + g r ^^^w''{lo's)oQ{s)ds. 
Jo s 
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Let us set the following initial conditions Q(0) = and Q{0) G S^j. We are interested in 
finding a solution on the interval < z < L. What we have here is a stochastic Volterra 
equation with (Fredholm) kernel 



s . 



(44) 



which is continous and 

Wk'^iz, s)\\h,-, < gM X[o,.](s)s-i 1^ _ ^1 ^ (45) 

Now we have to see what are the conditions that make equation (^3D solvable. We 
propose as ansatz the usual resolvent for convoluted kernels, that is, 



Kh{z,s) = Y,kP{z,s), 



(46) 



n=l 



such that 



Q{z) = Qqz+ / Kh{z,s) o [Qq-s] ds 



QqO 



z+ Kij(z,s)sds 



with the K^^^ given inductively by 

K^H^^\z,s)= I K^H\z,u)ok"{u,s) du, with n > 1, 

J s 

= k^'iz.s). 



K^i\z,s) 



(47) 

(48) 
(49) 



0ksendal|T^ found this is the unique solution for bounded kernels in the distribution 
Hida space. The same procedure can also be applied to fractional Hida spaces. But it 
can be seen from (^51) that it is not our case when s — > 0. Nevertheless, a new norm can 
be defined so that under it the former resolvent exists. Let it be. 



\k 



Hi 



\LP'P'{J),-g 



sup / / \\G{z)k" {z, s)F{s)\\H,-qdz ds 

,.n._„ <1 Jj Jj 



(50) 



\\LP{J)-q 



where ||F||LP(j)^_q = (/^ \\F{s)\ffj_^dsY^'' is the Lebesque integral, and J = (0, L]. This 
norm is simplified using the Holder inequality the following relation can be proved: 

i/p' 



\k"\\Lp,p'{j),-q < min 



J \Jj 



\k''{z,s)r^^_/s 



j^(^jjk^{z,s)\\i^_^dzy\s 



p'/p 



dz 

i/p 



(51) 



Gripenberg et al p3| discuss the deterministic counterpart of this construction. They 
proved that a resolvent solution exists whenever the norm of the kernel is less than 
one. This theorem can be tracked back to our norm in the stochastic case. Hence, the 
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same hypothesis apphes for this stochastic Fredholm kernel: \\k^\ 



LP'P'{J),-g 



< 1 for some 



g > 0. Then applying equation to the bounding condition (^) we find 



i/p( 



i/p' 



< 1, 



(52) 



\k"\\LV,v'(J),-q<gMp'' . 

since M is a small constant and g '^1. This guaratees the convergence of the proposed 
ansatz. 



The solution represented as a series of convoluted kernels, (|47|)-(^9D, is useless for 
calculations. Next, we will prove that a fractional chaos expansion exists for the solution. 
Let us take the second term in the Wick product of equation (^7|) , it can be written 



X(z)=z + 



n=l 



s ds 



n=l 



K^^\z^ s)s ds 



(53) 



because it converges absolutely. The general term in this series can be written. 



K^^\z, s)s ds = qIq 



l~H 



[S2 - Si)(si - 0) 



dso 



Sl 



Si 



oW"{si)dsi 




[S2 - Si)(Si - 0) 



'\z,S2) dsidB 



J sx 



z pz rz 



{gl 



1-H\2 



1-H\n 




J siJ S2 



[S3 - S2){S2 - Si){si - 0) 
S2S1 



dss dB^JB^ 



'S2 Sl 



Sn-1 



[z - Sn){Sn - Sn-l) ■ ■ ■ (gj - 0) _ _ _ 



SnS'n—l ' ' ' Sl 



[Z - Sr 



n 

i=l 



\Si - Si-i, 



■X[s,-^,z){.Si) 



H 



dB" ■ ■ ■ dB^^ 



zigz^ll-^'rf f^"\sn, . . . , 3i) dBZ ■ ■ ■ dBl 



(54) 



We have used property (|29| ) to build the above adimensional integrals, and defined 



n 

Sl) = (l-s„) JJ 

i=l 



/^"nSn,...,Sl) 



[Si — Si^i 



-X[s,_i,l)lSiJ 



(55) 



with So = 0. Now, the symmetrized form /'^"^(s„, . . . , Si) = (1/"^!) J2aeu f^"Ksa„, ■ ■ ■ ,Sc 
induces the following relation 



/ /(«)(s„, . . . , Sl) dBZ ■■■dBZ= [ /^"H^n, . . . , 3i) dBZ ■ ■ ■ dBZ, (56) 



and finally. 



x{z) = z\i + j2 \rf^''\sn,...,si) 



n=l ^ ^+ 



dB^ ■■■dB 



H 

Sl ( ' 



(57) 
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where g = IqqIz/Iq)^ . This will be nothing else but the fractional chaos expansion 
provided 



(n)|2 



< OO 



(5^ 



n=l 



holds. In fact this condition express nothing else that the existence of the variance of 
the process, 



(n)|2 



n=l 



(59) 



The search of an upper bound for the sucession of 0-norms, given that the are 
symmetric, is straightforward: 



1/ 



(n)|2 



/»(S„,...,S1)/»(4,,...,4)X 



X <P{Sn, SJ 



)(si, s[) dsn ■ ■ ■ dsi ds'^ ■ ■ ■ ds[ < 
- / /■■■/ / 4>{sn,s'J---(j){si,s[)dsn---dsids'^---ds[, (60) 

because of definition ( p5D and the fact < Sj — < Si (idem < — s[_^ < s[) the 
last inequality follows. Observing that 

"1 /•! /•! /•! 



/ / (j){sn,s'Jdsnds'^ = H{2H -1) / \sn- s'J^" dsnds'^ 



1 

2 L 



1 — Sn-l) + (1 ~ -^n-l) " ~ s[ 



,2H 



< 1, 



(61) 



we iteratively apply it in ( |60| ) to find: [/'•"•'l^ < 1. Thus, the chaos expansion exists for 
all z < L whenever 

H 



(62) 



is satisfied. 



4. Ray-light Statistics: An Example 

In this section we will use the stochastic ray-equation solution to study the statistical 
properties of the turbulent refractive index. Both coordinates of displacement are 
independent, and they also hold the same (non-coupled) differential equation. There is 
enough to consider the 1-dimensional then. The parameter election (p!6D , we have used 
in our treatment, also defines the meanining of the transversal velocities, for they are 
the angles of deviation. Being the velocities continuous we can set. 



(63) 



Since our solution is dependent of the initial refractive angle 6, its behavior at the 
boundary, e ^ 0, must be known. This boundary is just the interface between turbulent 
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and resting air. Henceforth, we will also model the initial angle as a fractional Brownian 
motion, 



c[ X[o,e)(s)rf5f = ci?^(e) 



(64) 



the constant c is adimensional and mesures the strength of the noise. The length e 
functions as a kind of correlation distance, as it goes to zero we are examining the 
properties of the interface's short-range correlation. 

Besides, any stochastic process can be put in terms of the spans described in the 
appendix, and these depend on the construction of stochastic integrals by step functions. 
So, even if the former model needs to be corrected-maybe the interface introduces long- 
range correlations-the next results are useful; since, they are the building blocks for 
more complex stochastic processes. 

The solution (|47|) is written using the chaos expansion (|57|) and the initial 
conditions: 



e{e)oX, 
zcB" 



Z 

\ n=l -^^l 



f^-\s^,...,s,)dBl---dB^^ 



H 



(65) 



From the Wick product properties is easy to see that 



E[Q{z)] =zcE[B 



Hi 



X E 



n=l -^^l 



= ■ E[l] = 0. (66) 

The evaluation of the variance from experimental data is the most common topic in 
many works related to the optical properties of turbulence because it is directly related 
to the Structure Constant. Hence, we calculate it using property (|A.39|) , 



E(D^ X,)2 + E(Xf)|x[o,e)| 



(67) 



E(X^) was already evaluated in the last section. The fractional Malliavin derivative 
appearing at the right-hand side demands elaboration, property ( [A.37| ) implies 



{DtX,)X[o,e){^)ds, 



the 0-differential is linear so we have 

oo 

Dtx, = zY,rD, 



n=l 



[ f^-\sn,...,s,)dBZ---dBf^ 



(68) 



(69) 



We are going to compute these derivatives now: let us fix n > 2, from the first theorem 
( |A.42| ) we can commute the stochastic integral and 0-differential, 



Dt 



f 



in). 



s,)dBZ---dBZ 



[ Dt [ f^-UBl-.-dBl 
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+ / /» dBZ ■ ■ ■ s^)ds,. (70) 

We can recursively commute the operators, and each time we do so another integral 
as the last in the right-hand side of the equation above is added. When we reach the 
innermost integral, we use property (|A.38|) to arrive to 



/ /»0(.,.„)c/5„c/i?2_^---ci< + 
[ f^^UBlcl>{s,Sn-i)dsn-i---dBl + ...+ f f^'^UBl---dBlc^{s,s,)ds, 



n 



/» ct>{s.Sn)dSn 



dBZ_,---dBi, 



(71) 



the last step is because of the symmetry of /*^"\ Instead, for n = 1 we just use property 

(EE)-- 



Dt 



f^'\si)dBl 



f^^\s,)<P{s,s,)ds,. 



Then we build the fractional Malliavin derivative from the series 



(72) 



D^^^^^^^ X,= z~g{ I f^'\s') X[o,e) (s)0(s, s') ds' ds + 



oo „ 
n=l -^^l 

its second moment is 
E 



H 



dBf ■ ■ ■ dB^^ 



2-2 

z g 



f^^\s') X[o,e)is)(f){s, s') ds' ds 



+ 



+ J2{n + 1)2^2" 



n=l 



X[O,e)(s)0(s,s') ds'ds 



(73) 



(74) 



This series converges, we apply the same procedure as before to find a bound for the 
integrals. What is more, each norm appearing in the series is bounded by the zero term. 



f^+'Hs',-) X[O,e)(s)0(s,s') ds'ds 



< 



/ (1 — s')0(s, s') ds' ds 
IJo Jo 



{2(2/1 + 1) 



e2^^+i _ (1 _ e)2^^+i] - 



H 



2(2/f+l)L^ ^ J 2/f+l^ ^ 2} - 4 ^^^"^ 

Thus, the existence of (|73D is guaranteed. Finally, the variance of the displacements, 
with the norm 



\X[o,e)\l = H{2H-l) 



\u — s\ duds = e , 



(76) 



^0 
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is written as, 



2/ \~\ 2 2 

'z] \ = Z C 




f^^\s') X[o,e)(s)0(s,s') ds'ds 



+ 



+ E 

n=l 



~2n 



f 



X[o,e){s)4>{s,s') ds'ds 



(77) 



Now, as the correlation distance goes to zero we recover the initial condition. While 
terms coming from the second moment of Xz banish (they are all bounded and multiplied 
by e^^), it is not the case with those coming from the fractional derivative. We will 
not go throught copious calculations since we are interested in a general outline of the 
solution; thereof, the solution can be expressed as 



z^c^g'^ 



+ 



+z\^f f\n + 1)^ 



ra=l 



/^"+'^(s',-) XW){s)4>{s.s') ds'ds 



(78) 



5. Remarks & Conclusions 



In this paper we provided a fractional Brownian motion model for the turbulent index of 
refraction; afterwards, we used this model to build a stochastic ray-equation-a stochastic 
Volterra equation. We have given a (unique) solution. Our analysis just covers the 
H > 1/2 case which is meant for non-isotropic or near diffusive range turbulence. 
These Hurst exponents disagree with those employed in optical experiments, for the 
turbulence is assumed nearly isotropic and homogeneous. We pretend to give a glimpse 
of the if < 1/2 problem from the solution we have found. 

In particular, this solution strongly depends on the initial conditions; our election 
of the initial angle at the example is the natural choice given the behavior of scalar 
quantities in turbulence. Under this condition, we apply solution ( [78D to estimate the 
centroid's variance of a laser beam at a distance L: 

E[Q{Lr] oc C'//'-'^L'^'^[l + rm, (79) 

and the series of (/)-norms is bounded, F{g) < g^/{l — g'^Y, whenever [L] ~ 1 m. 

Now, as if — > 1/2 the variance of the displacement approaches to C'^l^^^^L^. We 
notice that this behavior is the same found in Consortini a/ 0. But it not comes 
from the isotropic Kolmogorov's law, with the hyphoteses made at the introduction it 
corresponds to a Brownian motion {H = 1/2). That is, given a gaussian process with 
Structure Function's power as in (|1|) the result from Consortini does not follow. We have 
shown this can only be achieved on inert ial- diffusive conditions, for C„ \ 1/3; shortly, 
the Beckman's hyphoteses [^] do not apply to gaussian models used in Fluis Dynamics. 

Nevertheless, the approach we have introduced is just the begining of a long journey, 
since we must examine the stochastic ray-equation when H G [1/3, 1/2). By doing so. 
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we will confirm that the power law for the variance of the displacements is univocally 
determined and no discontinuities arise, equation (|79D is still valid. But, this will require 
the introduction of other tools-even new for the Stochastic Analysis-since in this range 
any smoothness property of the fBm processes is lost and this inquiry will be the 
topic of future works. 

A step further should be considered afterward. The proposed model for the 
turbulent index of refraction ( p6|) just aproximates the Structure Function. After we 
give a solution to the stochastic ray-equation for the whole range of Hurst paremeters 
with this model, we will start examining other functionals of the fractional Brownian 
motion, which refines our theoretical results against observed properties. 
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Appendix 

Here we will introduce briefly the elements needed to build a stochastic calculus for . 
The reader can find a complete reference in Hu and 0ksendal ||2^ and Duncan et al 

Let H G (1/2, 1) be a fixed constant, and let us define 



(l){s,z) = H{2H -l)\s~ zf^~\ s,zeR. (A.l) 
Then / G L|(]R), if it is measurable and 

\f\l--= I f f{s)f{zms,z)dsdz<cK>. (A.2) 

Also the inner product can be defined in L|(M); 

{f,g)4,:= [ f f{s)g{z)4>{s,z) dsdz, for all /, ^7 G (M) , (A.3) 



and it becomes a separable Hilbert space. 

Let iS(M) C L'^i^) be the Schwartz space of rapidly decreasing smooth functions on 
M. From its dual the probability space Q = iS'(M) can build, it is the space of tempered 
distributions u on M. Its associated probability measure, fi^, can be found applying the 
Bochner-Minlos theorem, 

E[e'<-'^>] := [ e*<-'^>rf/x^(cu) = e-^l^l^ (A.4) 
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where {u, f) is the usual pairing between elements in the dual and functions on M. With 
this definition it is easy to prove that 

EK-,/)]=0,andE[(-,/)2] =1/1^. (A.5) 

The triplet B{VL) , /i^) is thus a probability space-i3(f2) is the borel algebra 
on i7-usually called fractional white noise probability space. Then let LF'{ii^) = 
L'^{Q, B{Q) , fiff)) be the space of all the random variables X : ^ M such that 

||X||i.(^^):=E|X|^<oo. (A.6) 

Hence, those functions in define the set of random variables of the form 

f{u) = {u!,f) which is included in ^^(/i^); that is, the condition (|A.2|) induces square 
measurable random variables because of equations (|A.5| ). 

It can also be shown that S{M.) is dense in L^(M). This means that for any 
/ G we can build a series /„ G M such that /„ — > / in What is more, the 

following limit exists in L^(/i^): 

lim(o.,/„):=(cu,/). (A.7) 

n— »oo 

We can now define the fractional Browinian motion process as follows: 

B'^iz) := B'^iz, u) = {u, X[o,.)) G L\fi^). (A.8) 

In this definition B^ is thought to be the 2;-continous version of the rightmost hand side 
term. The step function X[o,z) '■ ^ 1] is defined by 

{1 iiO<s<z 
-1 iiz<s<0. (A.9) 
otherwise 

So taking / G L|(R), approximating it by step functions, and then using property (|A.7| ) 
traduces definition (|A.8| ) to 



{uj,f)= / f,dB''{z,u). (A.IO) 

It can be verified using the same procedure for /, (7 G that 

n{uj){uj,g)] = {f,g)^. (A.ll) 
Again let / be as above, as it is shown in Duncan et al once defined 

^(/)=exp(^^/rfi?^-^|/|j), (A.12) 

the linear span 

£=\Y, (^kSifk); n G N, afc G M, /fc G (M) for G {1, . . . , n} i , (A.13) 



, fc=i 



is dense in L^(/i0). 
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Also there is another functional expansion from where this Lebesgue space with 
measure fi^p can be build. This expansion is useful, in particular, to introduce some 
tools we will use through this work. To do so let us define the Hermite functions|p7|| 

U^) = {2^-\n-l)\V^)-'^\-^'/^H^_,ix); n = l,2,--- (A.14) 

on M. This set of functions is an orthonormal basis in L^(]R). We can map this basis to 
an orthonormal one = r^^.^„ G L^(M), by means of the isometry 



r^fis) = CH I {z- sf-'/'f{z) dz, (A.15) 
,00) 



where 



H(2H -1)T(^ - H) 

It can also be shown that 

Us)Hs,z)ds = CH I {z-s)''-^'^Us)ds, (A.17) 
J (—00,2;] 

because the ^„'s are an orthonormal basis these integrals are also smooth. 

Let X be the set of all finite multi-indices a = (ai, ■ ■ ■ , a^) of nonnegative integers, 
we define 

n^icu) :=H^,{{cu,ii))---H^J{uj,U)- (A.18) 
In particular, if we let := (0, ■ ■ ■ , 0, 1, 0, ■ ■ ■ , 0) denote the i'th unit vector then we get 
from (|A.10|) and the definition of Hermite polynomials 

HAu^) = m^,l)) = = [ Us) dB^. (A.19) 



These functionals are elements of L'^{fi^), and they form its basis. That is, for 
X G L'^{fi(i)) there are Cq, G M and a G X, such that 

X{uj) = Y,Cana{uj), (A.20) 

and also 

These coeficients are given by Cq, = E[X7Ya] /ctl. With this property at hand we now 
define the fractional Hida spaces: the fractional Hida test function space Sh as the set 
of all 

'?/'(a;) = aa7-^Q,(u;) G L'^{fj,(p), such that 

(A.22) 

MH,k = E ««(2 N)^° < 00, for all k G M, 

where 

(2 ny = Yli2jy^ for any element 7 = (71, ■ ■ ■ , 7^) G X. (A.23) 
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The fractional Hida distribution space S*^ is the set of all formal expansions 
Y{ijo) =^^6/37-^/3(0;), such that 

(A.24) 

W\\H,-q = ^(^^- «J(2 N)"^^ < cx), for some q en. 

With these definitions is not hard to see that Sh C L'^{fi(f,) C S^. 

It is now time to show how the fractional white noise and integration with 
respect to can be defined. Let us first calculate the expansion for the stochastic 
integral in ( [A.10| ). For any / G L^(]R)-any given deterministic function- we have from 
equations ( pOOD and ( |7rTl| ): 
„ 00 

/ fsdB^ = Y.{fAi)^nA^)- (A.25) 
When / = X[o,z) in the left hand side we recover (|A.8|) and the following relation holds 



ik{s)4>{s, u) ds \ du 



(A.26) 



if we check its norm 



E 

k=l 



'Cfe('5)0(s, u) ds J du 



[0,2) \JR 
00 



k=l ^ ^ 



(A.27) 



(C is the Riemann's zeta function) because 

4(s)0(s,m) ds = (m - s 

J { — OD,u] 



.H-3/2 



^k{s)ds 



(A.28) 



according to Szego |28|. When g > 4/3 the former inequality shows that B^ is continuos 
and differentiable in S*fj. Hence, 

j-^B^{z) = (^jj,{s)cP{s,z)ds^ HA^) ■■= W'iz) e S*H (A.29) 

is the formal definition of fractional white noise. Again it is also continous, when z > s 



\W«(z) - W«(s)\\l_, 



k=l 

00 



^k{u)(f){z,u)du- / ^k{u)(f){s,u)du 



k=l 



[^,^s)-u]''-'/'Ms + u)\du 



{2k)-'' < 



{2k)-'' < 



<2-«4M2C(g + - 



6 



[{z -s)-u] 



H-3/2 



22-'?4M2 



(2if- 1)2 



C q 



1 



iz — s 



,2H-1 



(A.30) 
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We now define the Wick product, like in Hu and 0ksendal, as follows: let 



a,l3eJ 7GX \o+/3=7 / 



{XoY){u) = 

this product is commutative, associative and distributive like the usual product in M. 
Under the norm || • Wn-q it can be shown that all the power series defined in the real 
line have their counterpart in the distribution Hida space. If we understand by X*" the 
n-times Wick product of X G then the Wick exponential is 

00 ^ 

exp^(X) = J2 (A-32) 

n=0 

In particular when X = {u, f) = fsdB^ holds, then 

expO((c.,/)) = ^(/). (A.33) 

The right-hand side was defined in (|A.12| ) , and therefore we have shown the link between 
the two representations presented up to now. 

It is appropiate to introduce the fractional Malliavin derivative or 0-derivative, for 
X e L\^i^) and G ^^(M), 

D^gX{uj) = \iui-\x{uj + 5 [ {^g){u)du))~ X{uj)\ (A.34) 

where {^g){z) = J^(f){z,u)gudu. Furthermore, if there exists a function DfX such that 

D^gX = [ {DfX) gs ds, \fg G ^^(R), (A.35) 

we say that X is 0-differenciable, and DfX is the 0-differential. These are differential 

operators, and they also present the following properties: let X be as always and 
f,g:R^R then 

D^JiX) = f'{X)D^gX, (A.36) 

D^g{ujJ) = {f,g)^, (A.37) 

Df{u;,f)= [ 4>{u,s)fudu (A.38) 



Moreover, it can be proved for any two X, F G L^(/i(^)-decomposed by the span 
using the Wick product and the (/)-derivative properties that. 



E 



Xo / 



Yo [ g^dBf 
Jm. 



E[iD^fX)iD^gY)+XYif,g)^] . (A.39) 



These properties allows to change the integrator inside ([A.10| ) by a stochastic function 
X : M X f2 — > M. The basic procedure consists of building a Riemann sum, replacing the 
standart product by the Wick one, 

n-l 

^„(X) := J]X,^ o (5^^^ - (A.40) 
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and finally observing its convergence (whenever E|X|^ < oo) under the E|-|^-norm. 
Moreover, the following equality holds 



X, dB^ = XsO ds, (A.41) 

[0,L) J[^,L) 

while the integral on the left-hand side represents the limit of ( |A.40| ), the right-hand 
side is just the integral evaluated under the Hida expansion of the Wick product defined 
in 

We will complete this appendix enumerating some useful theorems from Duncan et 

al 

Let Xz a stochastic process defined as above, anc? sup^g^,/.) | -^2^ I ^ < Also, 
let r]z = fs dB^ . Then for s, 2; G [0, L) 

Dtr]z= [ DfX^dB^+f X^(f){s,u)du. (A.42) 

J[0,z) J[0,z) 

Let Tjz = /jg fs dB^ and F{z, x) : ]R_|_ x M ^ M, where f is continous and F has 
second continous derivatives. Then 



f dF f dF 

F{z,7]z) = F{0,0)+ (s,T]s)ds+ {s,T]z)f,dB^ 

J[Q,z) OS i[o,^) dx 

+ / TTY^^^^^^ / <P{s,s')fs' ds'ds. 

^[0,2) ox J,Q ) 



'[0,2) '-'•^ J[0,s) 

If X E L'^{ji^) then there exists a secuence {/„ G L'^^{W\_)}nm such that 
T.n=i \fn\l < 00 and 

00 „ 

X = E[X] + V / /„(si, ■ ■ ■ , s„) dBf^ ■ ■ ■ dBl. (A.44) 

n=l 

Here it is understood L'^^{W\) as the n-dimensional space of symmetric functions and 

\fn\l = / fn{Si, Sn)fn{s[, " " " , <) 0(Si, s[) ■ ■ ■ 0(s„, s'J dSi ■ ■ ■ dSn ds[ ■ ■ ■ ds'^. (A.45) 
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